Validity of Effective Potentials in Crowded Solutions of Linear and Ring Polymers with Reversible Bonds

We perform simulations to compute the effective potential between the centers-of-mass of two polymers with reversible bonds. We investigate the influence of the topology on the potential by employing linear and ring backbones for the precursor (unbonded) polymer, finding that it leads to qualitatively different effective potentials. In the linear and ring cases the potentials can be described by Gaussians and generalized exponentials, respectively. The interactions are more repulsive for the ring topology, in analogy with known results in the absence of bonding. We also investigate the effect of the specific sequence of the reactive groups along the backbone (periodic or with different degrees of randomness), establishing that it has a significant impact on the effective potentials. When the reactive sites of both polymers are chemically orthogonal so that only intramolecular bonds are possible, the interactions become more repulsive the closer to periodic the sequence is. The opposite effect is found if both polymers have the same types of reactive sites and intermolecular bonds can be formed. We test the validity of the effective potentials in solution, in a broad range of concentrations from high dilution to far above the overlap concentration. For this purpose, we compare simulations of the effective fluid and test particle route calculations with simulations of the real all-monomer system. Very good agreement is found for the reversible linear polymers, indicating that unlike in their nonbonding counterparts many-body effects are minor even far above the overlap concentration. The agreement for the reversible rings is less satisfactory, and at high concentration the real system does not show the clustering behavior predicted by the effective potential. Results similar to the former ones are found for the partial self-correlations in ring/linear mixtures. Finally, we investigate the possibility of creating, at high concentrations, a gel of two interpenetrated reversible networks. For this purpose we simulate a 50/50 two-component mixture of reversible polymers with orthogonal chemistry for the reactive sites, so that intermolecular bonds are only formed between polymers of the same component. As predicted by both the theoretical phase diagram and the simulations of the effective fluid, the two networks in the all-monomer mixture do not interpenetrate, and phase separation (demixing) is observed instead.


INTRODUCTION
Single-chain nanoparticles (SCNPs) are soft nano-objects, of size in the range 3−30 nm, which are synthesized through purely intramolecular cross-linking of functionalized polymers (precursors). 1 Both for their size and for their internal malleability that allows for quick response to environmental changes, SCNPs are promising macromolecules for applications as catalytic nanoreactors, drug delivery nanocarriers, and biosensing probes, to name a few. 2−8 Depending on several factors implemented on the precursor, such as the solvent conditions, its molecular topology, chain stiffness, or the presence of crowders, the resulting SCNPs present a broad range of structural conformations, from very sparse objects 9,10 to more compact and even nanogel-like SCNPs. 11 −13 In the usual route (linear precursors in good solvent), the resulting SCNPs are sparse objects where short-range loops dominate the distribution of cross-links. This is a direct consequence of the self-avoiding conformations of the linear precursor in good solvent. In such conformations contacts between monomers separated by long contour distances and formation of long-range loopswhich are efficient for folding into globular shapesare infrequent. The SCNP conformations are dominated by short loops and have scaling exponents of ν ∼ 0.5 for the dependence of the size on the number of monomers (R ∼ N ν ), far from the globular state ν ∼ 1/3.
Synthesis of SCNPs has been traditionally dominated by irreversible intrachain cross-linking of the precursor. In recent years, growing efforts have been dedicated to broaden the functionalities and areas of applicability of SCNPs through the implementation of reversible bonds in their backbone via noncovalent and dynamically covalent interactions. The current SCNP chemistry toolbox of reversible bonds includes, among others, hydrogen bonds, metal complex formation, hydrazone, enamine, anthrazene, and so on. 14−17 Because breaking and formation of these bonds can be activated and tuned through factors such as temperature, pH, or light, the single-chain character of these objects in solution is lost when their concentration is high enough, leading to the formation of aggregates and eventually a percolating network. Because of the reversibility of the bonds the bonding pattern of the network is dynamic, allowing for viscous flow of the material and for physical gelation if the external stimuli are switched off (e.g., by decreasing temperature). The possibility of designing smart polymers that can reversibly transform from solutions of SCNPs to hydrogels has been demonstrated experimentally. 18,19 These findings pave the way to use polypeptidebased SCNPs as building blocks for biocompatible and biodegradable materials with self-healing properties and applications in tissue engineering. 19 Recent simulations have investigated the transition from a solution of sparse SCNPs at high dilution to a dynamic network in semidilute and concentrated conditions for a system of linear chains with reversibly bonding sites in their backbones. 20 Some remarkable results have been reported: (i) the intramolecular bonds still form the majority of the overall bonding, and the connectivity of the network is mediated by a few intermolecular bonds per chain; (ii) the bonding pattern of the network is dynamic, and the polymers can diffuse long distances through breaking and formation of bonds at different sites without losing their connection to the percolating cluster; (iii) the size and shape of the SCNP conformations at high dilution are essentially unaffected by crowding and remain in the network even at densities far above the overlap concentration. The latter is a rather unusual result, in clear contrast with the shrinkage found for other sparse objects such as simple (unreactive) linear chains and rings, which by increasing the concentration change their conformations from self-avoiding (Flory exponent ν ≈ 0.59) to random walks (ν = 1/2) in the case of linear chains 21 and to fractal ("crumpled") globules (ν = 1/3) in the case of rings. 22,23 The weak effect of crowding on the molecular conformations of the reversible SCNPs is inherently related to the formation of intermolecular bonds. Indeed, when the SCNPs are prepared at high dilution through irreversible intramolecular cross-linking and are transferred to high concentration, with no intermolecular bonding, they show a collapse similar to that of rings to crumpled globular conformations. 24, 25 The effective potential between two macromolecules separated by a given distance is the free energy needed to bring them from the infinity to that distance. Unlike in hardcore colloids, the free energy cost for full interpenetration of the macromolecules (zero distance) is finite because their centers-of-mass can coincide in space. The cost of full interpenetration strongly depends on the topology and internal deformability of the two macromolecules, typically varying between a few and tens of times the thermal energy. 26,27 Averaging out the molecular internal degrees of freedom and keeping one or a few relevant coordinates (usually the centersof-mass) reduces the system to an effective fluid of ultrasoft particles interacting through the effective potential. 28−32 This methodology allows not only for simulating much larger scales than in the monomer-resolved models but also for the treatment of the system by methods from liquid state theory, producing a powerful tool for predicting large-scale organization and phase behavior. 33,34 A major limitation of this approach is that because the effective potential is derived for a pair of polymers in the absence of others, it neglects the manybody interactions that are present in a crowded solution or a melt. This approximation works well below and even slightly above the overlap concentration (i.e., the concentration at which the mean intermolecular distance is of the order of the unperturbed molecular size). However, it fails dramatically far above the overlap concentration when many-body effects become a dominant contribution (shrinkage of molecular size being a manifestation of them). A paradigmatic example is the nonemergence of the cluster crystal phase predicted by the effective potential for flexible ring polymers, 35 which instead collapse to crumpled globular conformations that hinder the full interpenetration required to form the cluster nodes.
As mentioned above, when linear polymers with reversible bonds assemble into a dynamic percolating network, they essentially maintain the SCNP conformations adopted at high dilution. 20 This result suggests that many-body effects can be negligible for this system, and the interaction of a tagged pair with their neighboring molecules is effectively given by a flat energy landscape not affecting the effective mutual force between the two polymers of the tagged pair. In such a case, the validity of the effective potential to describe the static correlations between molecular centers-of-mass could extend to unprecedented densities far above the overlap concentration. With this idea in mind we investigate, by means of simulations, the validity of the effective potential for a system of generic bead−spring polymers that switches from a solution of SCNPs at high dilution to a reversibly cross-linked polymer network at high concentrations. We explore a broad concentration range between both limits as well as the effect of the molecular topology of the unbonded polymer (linear or ring) and the sequence of reactive sites (with different degrees of randomness) along the molecular backbone. We test the accuracy of the effective potentials by comparing simulations of the real all-monomer systems with their corresponding effective fluids of ultrasoft particles as well as with predictions from the test particle route. 36 We also test the approach for a ring−linear mixture as well as for a two-component linear/ linear mixture with orthogonal bonding chemistry, where intermolecular bonding is only allowed between chains of the same component. We find that both the topology of the unbonded polymer and the specific sequence of the reactive sites along the polymer backbone have a strong impact on the effective potential. As suggested by the weak effect of the concentration on the size and shape of the linear polymers with reversible bonds, the simulations confirm that the effective fluid provides a very good description of the real system at densities far above the overlap concentration. In a similar fashion to the case of ring polymers with no bonding, the effective fluid approach is less satisfactory for the ring-based system, and the predicted clustering behavior is not found in the real system. The effective potential becomes much more repulsive when intermolecular bonding is switched off. As a consequence, the effective binary fluid representing the mixture with orthogonal bonding chemistry shows demixing. This behavior is confirmed in the all-monomer real mixture, which shows spontaneous demixing within the simulation time scale. This striking result suggests that experimental interpenetrated networks with reversible bonds are kinetically trapped states where demixing is prevented by large barriers arising from long bond lifetimes and entanglements.
The article is organized as follows. In section 2 we define the model and interactions implemented in the all-monomer simulations. We also give the simulation details for the computation of the effective potentials and briefly describe the analytical test particle route approach. In section 3 we report a critical analysis of the obtained effective interactions as a function of the topology of the precursor and the specific sequence of reactive sites. In section 4 we present theory and simulation results for the solutions at different concentrations and for the phase behavior of the mixtures and discuss the validity of the effective potentials to describe the behavior of the all-monomer real systems. In section 5 we summarize our conclusions.

MODEL AND SIMULATION DETAILS
The precursors are modeled as fully flexible linear chains or rings made of 200 beads (monomers). A fraction of these monomers f = N r /N m = 20/200 = 0.1 are reactive, where N r and N m are respectively the number of reactive sites and the total number of monomers. The reactive sites can form and break bonds with other reactive sites within the same polymer (intrabonding) or with reactive sites belonging to other polymers (interbonding). In all cases, we perform Langevin dynamics simulations. In the first step we use them to obtain the effective potentials (section 2.2); subsequently, we use them to simulate the effective fluid at different concentrations, and we compare results with simulations of the corresponding all-monomer system (section 2.3). Moreover, we compare simulation results with theoretical calculations by the test particle route (section 2.4).
2.1. Model. We describe the polymer chains by the bead− spring model of Kremer and Grest. 37 Thus, excluded volume interactions among the beads are modeled by the Weeks− Chandler−Andersen (WCA) potential: The permanent bonds leading to the connectivity of the precursor are implemented via a finite extensible nonlinear elastic (FENE) potential between consecutive monomers. This is given by where K F = 15 and R 0 σ = 1.5σ is the maximum elongation of the bond. For implementing the reversible bonds between the reactive sites, we adopt the potential introduced by Rovigatti In our system we set the capture radius r c = 1.3σ while ϵ ss = 12ϵ and K ss = [σ/2(σ − r c )] 2 . With these choices the potential and force are continuous and zero at the capture radius. Moreover, the potential is short-ranged and has a deep attractive minimum of energy V ss min = −12k B T (which can be seen as the bond energy) at the distance r min = 1.0σ. When the distance between two reactive sites is smaller than r c , the interaction of eq 3 becomes nonzero and attractive and the sites form a mutual bond. The bond is broken when a fluctuation moves the mutual distance beyond r c . Because we wish to limit the valence to a single reversible bond per reactive site, we make use of the swapping algorithm introduced by Sciortino. 39 Thus, we add a repulsive three-body contribution in such a way that it is switched on when a reactive site k enters the capture radius of a reactive site i that is already bonded to another one j. The three-body potential is defined as where the sum includes all i, j, k triplets, and Therefore, 0 < V 3body (r) ≤ |V ss min | for each triplet, and when a triplet is formed, the energy decrease resulting from the new bond is compensated by the three-body repulsive term without changing the potential energy of the system. This three-body term makes triplets very short-lived and spontaneously leads to bond swapping, which speeds up the exploration of different patterns of the bond network. Moreover, monovalent bonding is governed by a Hamiltonian, unlike methods based on random choices when one site can bind to more than one candidate. 9 Because the polymer is fully flexible, a monomer qualitatively corresponds to a Kuhn length (of the order of 10 monomers), so that the actual fraction of reactive sites qualitatively corresponds to 1% in real polymers. Moreover, because bonding is nondirectional (unlike in patchy models), a reactive site qualitatively represents a functionalized pendant group with high flexibility. The former conditions are indeed common in experimental SCNPs, which are the natural state of our systems at high dilution. For simplicity, we set m = 1 for all monomers, so that the center-of-mass coincides with the Figure 1. Representation of the different interaction potentials used in this study. The combination of the WCA and FENE potentials results in a deep potential well that sets the mean length of the permanent bonds at r min,irrev ≈ 0.96σ. The combination of the WCA and the reversible bonding potential V ss (r) defines the mean length of the reversible bonds at r min,rev ≈ 1.0σ. The function |V ss min |Φ 3 (r) represents the contribution of a bond belonging to a triplet to the three-body potential (see text). geometrical center. The simulations were performed at temperature T 1.0 by using a Langevin thermostat with a friction coefficient γ = 0.05. 40 Equations of motion were integrated within the scheme of ref 41 by using a time step δt = 0.005. 2.2. Computation of the Effective Potential. The effective potential acting between the two polymers (1, 2) can be calculated by integration of the net force over the axis joining their centers-of-mass (see e.g. ref 29): F eff,12 = −∇ R 12 V eff (R 12 ), where R 12 is the distance between the two centers-of-mass. The net force experienced by one of the polymers is computed as the total force (nonbonded and bonded) exerted on its monomers by the monomers of the other polymer. In the following expression we consider the force exerted by polymer 2 on polymer 1: 12 , 1 (1), (2) where F i(1), j(2) is the force exerted by the jth monomer of polymer 2 on the ith monomer of polymer 1, the sum runs over all i(1), j(2) pairs, and the subscript on the right-hand side means that the average must be evaluated at the fixed separation R 12 . Obviously the expression accounting for the interactions of polymer 1 on polymer 2 just produces the opposite result, and integration leads to the same effective potential. The statistical averages of the components perpendicular to the axis joining the centers-for-mass are zero. We performed Langevin dynamics simulations where the positions of the centers-of-mass of the two polymers, and therefore their mutual distance R 12 , were kept fixed at every time step. We performed the simulation runs at the fixed distances R 12 /σ = 0, 1, 2, ..., 34, 35. For each distance we performed an equilibration run of 10 7 time steps, followed by a production run of at least 4 × 10 8 steps. To improve statistics as much as possible, the total force F eff was obtained by on-thefly averaging the summation of eq 6 over all the time steps of the production run. None of the initial bonds survived after typically 6 × 10 6 steps. Therefore, the simulations were long enough to achieve a good sampling of the ensemble of bonding patterns.
To test if there is any dependence of the effective interactions on the specific sequence of reactive sites along the backbone of the precursor, we consider three different cases to simulate for a couple of polymers with reversible bonds: (i) A random sequence of the 20 reactive sites with the constraint that there is at least n min = 1 nonreactive sites between consecutive reactive sites (to prevent trivial bonding). This case will be denoted as "gap1". (ii) A random sequence with the constraint n min = 4. This case will be denoted as "gap4". (iii) A periodic sequence; i.e., there is a constant separation n min = 9 between consecutive reactive sites. This case will be denoted as "periodic". In both cases i and ii the sequences of the two polymers are different, with the only condition that they have the same n min . Moreover, to assess whether even by using the same value of n min the specific realization of the sequences affects to the effective potential, we simulated two different couples (denoted as couple 1 and couple 2) for each of the cases "gap 1" and "gap 4". Figures S1 (linear) and S2 (rings) in the Supporting Information show the specific sequences of the simulated couples. Moreover, we investigated a couple formed by a linear chain and a ring. In this case the simulations were limited to the case n min = 1 (gap 1), and we used the first polymer of their corresponding couple 1.

Simulations of All-Monomer and Effective Fluids.
We performed Langevin dynamics simulations of solutions of the real all-monomer polymers and of the corresponding effective fluids. We explored the validity of the effective fluid approach in a broad concentration range below and above the overlap concentration, 21 which we define as ρ* = N m /(2R g0 ) 3 , where R g0 is the radius of gyration of the isolated polymer (i.e., in the absence of all intermolecular interactions). Therefore, if ρ = N p N m /V is the absolute density (number of monomers per volume), with V the volume of the simulation box and N p the number of polymers in the box, the reduced density (normalized by ρ*) is ρ/ρ* = N p (2R g0 ) 3 /V. In the case of a binary mixture of components (1, 2) we define the reduced concentration as ρ/ρ* = V −1 [N p,1 (2R g0,1 ) 3 + N p,2 (2R g0,2 ) 3 ]. Therefore, the overlap concentration of the binary mixture is For the isolated linear and ring polymers we find R g0 /σ = 9.93 and 7.72, respectively. Therefore, the density of monomers at the overlap concentration is ρ*σ −3 = 0.025 and 0.064 for the pure solutions of linear chains and rings with reversible bonds, respectively. For a 50/50 linear/ring mixture the overlap concentration is ρ*σ −3 = 0.036.
In the all-monomer simulations we investigated the pure systems of linear chains and rings with reversible bonds, a 50/ 50 linear/ring mixture, and a mixture of linear chains with orthogonal bonding. In the latter, bonding (intra-or intermolecular) was only permitted between polymers of the same component, and all WCA, FENE, and reversible bonding interactions were the same as in the other simulated systems, with the only condition that A-and B-reactive sites could not form mutual bonds and only interacted through the WCA potential. Although the breaking and formation of bonds can lead to concatenation of reversible loops in both the linear chain and ring-based systems, in the latter intermolecular concatenation between the permanent ring backbones must be avoided. Thus, nonconcatenated dilute ring-based systems were initially prepared and compressed to the target concentrations where they were further equilibrated. To prepare the linear−linear mixture with orthogonal bonding chemistry, configurations were taken from the one-component system and half of the chains were randomly assigned to each component of the mixture, which was further equilibrated with no intermolecular bonding between different components. Therefore, the final demixed state that we anticipated in the Introduction was reached spontaneously from an initially mixed state, demonstrating the robustness of this result.
The duration of the equilibration and production runs was typically 10 7 and 8 × 10 7 time steps, respectively. To improve statistics, eight independent runs were simulated at each concentration. The polymers moved at least 5 times their own diameter of gyration at all the investigated densities, thus guaranteeing a good sampling of the configuration space. In all cases the total number of polymers in the simulation box was N p = 108, with N m = 200 monomers and N r = 20 reactive sites per polymer, and the concentration was tuned by varying the box size. All the polymers had different random sequences of reactive sites corresponding to the case "gap 1". The effective fluids were simulated by using the corresponding effective potentials obtained for the couple 1 of the case "gap 1". Because of the much smaller number of degrees of freedom, in the effective fluids we simulated larger boxes of N p = 1000 effective particles, rescaling the box size to have the same concentrations as in the respective all-monomer systems. Tables S1 and S2 show the simulated box sizes for each allmonomer and effective system and the respective concentrations in absolute and reduced units. To test the effect of the box size, some concentrations in the effective fluid were also simulated with N p = 108 particles. Structural properties were not changed within statistics. This is demonstrated in Figure  S3, which shows representative results of the radial distribution function g(r) of the effective fluid of linear chains with reversible bonds. Data are shown at the lowest and highest investigated concentrations and in both cases for N p = 108 and 1000 effective particles (with the respective rescaling of the box size to produce the same concentration). No differences are found within statistics in the respective g(r)'s. Therefore, we conclude that finite size effects are not significant (except for the phase separating system of chains with orthogonal bonds, where the phase growth is obviously limited by the box size).
2.4. Test Particle Route to Fluid Structure. The test particle route (TPR) will allow us to compute the radial distribution function of the effective fluid by using the formalism of mean-field density functional theory (DFT) for inhomogeneous fluids. 36 In our study each particle of the effective fluid represents the center-of-mass of one polymer and interacts with the others via the effective potential V eff (r) computed as described in section 2.2. Within TPR, a particle is fixed at the origin of the system. As a consequence, the particle perturbs the system, and the density of particles around it changes from a constant bulk value ρ b to a spatially varying local density ρ(r). The external potential acting on the particle at the origin is equal to the effective potential, implying that the radial distribution function can be calculated as g(r) = ρ(r)/ρ b . Following the derivation from TPR based on DFT suitable for soft potentials (see the Supporting Information for details), the partial radial distribution function for the i, j components of a mixture of n components is given by 42 where ρ b,i is the macroscopic density of the i component, h ij (r) = g ij (r) − 1 is the ij component of the total correlation function, V eff,kj (r) is the interaction potential between species k and j, and the symbol * denotes convolution:

MOLECULAR CONFORMATIONS AND COMPUTATION OF THE EFFECTIVE POTENTIALS
3.1. Conformations of Two Interpenetrated Polymers. We have investigated effective interactions between two polymers with reversible bonds, namely two linear chains ("linear−linear"), two rings ("ring−ring"), and a linear chain and a ring ("linear−ring"). In all cases we have simulated two possibilities of bonding. In the first one (denoted as "all bonds") we carry out standard runs where the two polymers can form both intra-and intermolecular bonds. In the second case (denoted as "only intra") only intramolecular bonds are allowed; i.e., reactive sites belonging to different polymers only interact through the WCA potential and cannot form mutual bonds. Before discussing the effective interactions, we characterize conformations of the two interacting polymers through their radius of gyration R g and the asphericity parameter a. This parameter (0 ≤ a ≤ 1) measures deviations from spherosymmetrical conformations (a = 0) and is defined as where λ 1 ≥ λ 2 ≥ λ 3 are the eigenvalues of the gyration tensor of the polymer. Figure 2 shows for each of the topologies (linear, ring) and sequences of reactive sites (couples 1, 2 of gap1 and gap4, and periodic) the distributions of instantaneous values of R g and a collected from the trajectories. The data are shown for isolated polymers (mimicking the case V eff (r → ∞) = 0). Only the distributions for the first polymer of each couple of Figures S1 and S2 are shown. Figures S4 and S5 compare for each case the distributions of the two polymers of the couple. As can be seen, the ring polymers with reversible bonds are smaller and closer to spherical than their linear counterparts. For the same value of n min the specific sequences (4 in total for couple 1 or couple 2) have at most a minor effect on the distributions P(R g ) and P(a). However, Figure 2 shows that changing the typical distance between consecutive reactive groups does have a systematic effect on P(R g ). Namely, increasing n min leads to smaller sizes of the polymers. This is not surprising because longer distances between consecutive reactive groups promote the formation of longer loops, resulting in a stronger reduction of the molecular size with respect to the linear precursor. No significant effect of n min on P(a) is found. Figures S6 and S7 show the effect of the intermolecular interactions on the size and shape of the two polymers. The distributions P(R g ) (S6) and P(a) (S7) now correspond to a distance between centers-of-mass r = 3σ and allowing for intermolecular bonding. Similar results are found for other close distances. As can be seen in Figure S6, the mutual interaction tends to swell both polymers (the maxima of P(R g ) are shifted by about 15%) with respect to the isolated (r → ∞) case. The mutual interaction also tends to increase the asphericity ( Figure S7). A remarkable effect for the ring−ring case is that the two polymers do not swell in the same way. This can be explained by the fact that at short intermolecular distances one of the rings is threaded by the other one. Figure  3 shows typical conformations of the two polymers at mutual distance r = 0 (from (a) to (c): linear−linear, ring−ring, and linear−ring). Panels b and c illustrate the threading of one ring by the other polymer (this also occurs in the linear−ring case). The asymmetry found in the radii of gyration of the two interpenetrated rings is also reflected in their different asphericities ( Figure S7), though the effect is less pronounced than in P(R g ). Figure S8 shows the time dependence of the ratio of the instantaneous R g 's of the two polymers at a mutual distance r = 3σ. Orange curves are the bare data. Blue curves are the data smoothed by 100 point averaging. Whereas in the linear case the ratio quickly fluctuates, in the ring case it is relatively persistent, as expected for a threading mechanism. Moreover, the fact that the ratio for the two rings fluctatuates above and below 1 shows that both rings alternate their threading/threaded character, which is a signature of good configurational sampling. Figure 4 shows the effect of switching on and off the intermolecular bonds on the conformations of the two polymers at a close distance r = 3σ. In the case of the linear chains there is a tiny shrinking of the size of both polymers when intermolecular bonding is allowed, which is presumably due to the slight reduction of the fluctuations when a few intermolecular bonds connect the two polymers. A different behavior is observed in the pair of rings, whose sizes change in opposite directions when they form intermolecular bonds and their size disparity is reduced. Thus, the larger threaded ring  shrinks and the smaller threading ring swells. The combination of both effects, occurring in the asymmetric pair created by threading, reduces distances between segments of different polymers and facilitates the formation of intermolecular bonds.
3.2. Effective Potentials. In Figure 5a we show the effective potentials obtained for the interaction between two linear chains with reversible bonds. Figure 5b shows the corresponding results for two rings, and Figure 5c compares results of the former cases with the effective potential between a linear chain and a ring. All data sets in Figure 5c correspond to sequences gap 1, namely, the couples 1 of Figure 5a,b for the linear−linear and ring−ring case. For the linear−ring case the simulations used the first polymer of the couples 1 of the linear−linear and ring−ring cases. The symbols in all panels are the results obtained from the simulations. The solid lines are fits to a main function plus a tail, both given by generalized exponentials, 43 βV eff (r) = a 1 exp(−b 1 r m 1 ) + a 2 exp(−b 2 r m 2 ). The tail is added to obtain the best possible description of the data sets not only for the core of the potential but also for all distances and down to energies much lower than k B T. In general, the interactions between the linear chains with reversible bonds can be described by Gaussian functions (even without needing the tail), whereas exponents m i > 2 are needed for ring−ring and ring−linear interactions. Table 1 shows the functions that fit the potentials found for the linear− linear (all bonds and only intra), ring−ring (all bonds), and linear−ring (all bonds) interactions, namely, in the cases "gap 1, couple 1". These are the potentials that will be used in the simulations of the effective fluids discussed in the next section. Figure 5 reveals several trends. The potentials (both with intermolecular bonding switched on and off) are more repulsive for ring−ring than for linear−linear interactions, the linear−ring case being intermediate between the former two. This is consistent with the findings in the linear and ring precursors (i.e., in total the absence of both intra-and intermolecular bonding) and reveals that the topological interaction is again relevant. As can be seen in Figure 5a,b, if only intramolecular bonding is allowed (empty symbols), the amplitudes of the potentials are systematically higher than in their respective precursors (about 2.5k B T and 6k B T for linear− linear and ring−ring precursors, respectively 35 ). This result is not surprising because the presence of intramolecular loops, even if they are transient, enhances steric hindrance and topological constraints and creates higher effective barriers for interpenetration than in the respective precursors.
As can be seen, the effective potential becomes systematically stronger, with variations of about 30−50% in its amplitude, by moving from the "gap 1" to the periodic sequence of the reactive groups. As mentioned before, increasing the distance between consecutive reactive groups promotes the formation of longer intramolecular loops and reduces the molecular size. This hinders interpenetration and leads to stronger effective repulsions. For a f ixed value of n min the specific sequence of reactive groups has some small, but visible, effects on the effective potential (see e.g. data for the two couples "gap 4" in Figure 5a). We assert that this small effect will vanish for long polymers because pairs of segments interacting intra-or intermolecularly will sample a huge amount of local sequences within the scale of an interacting segment, so that averaging over the local sequences for a fixed n min will always lead to the same effective potential, irrespective  of the specific realization of the full sequence. On the other hand, we expect that the dependence on n min will persist for long polymers because n min affects to the typical intra-and intermolecular distances between reactive groups (e.g., a larger n min leads to longer intramolecular loops on average, which tend to increase steric repulsion). Figure 5 shows that when intermolecular bonding is switched on (filled symbols), the effective potentials experience a marked reduction with respect to the case of pure intramolecular bonding. Interestingly, the effect of the sequence of reactive sites when intermolecular bonds are allowed is the opposite to that found when they are not: increasing the distance between consecutive reactive sites decreases the effective interaction. As a consequence, the periodic sequences of reactive groups lead to the strongest reductions of the effective potential when intermolecular bonding is switched on (with differences of ΔβV eff (r = 0) ∼ −3 and −5 for linear−linear and ring−ring interactions, respectively).
The analysis of the number of bonds can shed some light on the origin of the former trends for the effective potentials. The average total number of bonds (intra-and intermolecular) is n tot ≈ 17 in all systems, i.e., 85% of the maximum n tot = 20 that would correspond to the fully bonded state. No differences are found within statistics, and this observation is independent of the topology of the precursor, the sequence of reactive sites, the distance between the centers-of-mass, and intermolecular bonding being switched on or off. Although the total number of bonds is unaffected, varying the former parameters leads to a different balance between intra-and intermolecular bonds. Figure 6 shows, for the cases of Figure 5 (same symbol codes), the variation of the number of intermolecular bonds, n inter , with the distance between the centers-of-mass of the two polymers. The number of intermolecular bonds increases by moving from gap 1 to periodic sequences, i.e., by increasing the distance between consecutive reactive sites. Because increasing such a distance eliminates the shortest intramolecular loops, the observed conservation of the average total number of bonds is achieved by exchanging the shortest loops by longer ones or by forming more intermolecular bonds. The second option is preferred, as shown by Figure 6. Figure S9 shows the distribution of instantaneous values of n inter at distance r 1,2 = 3σ. As can be seen, n inter can fluctuate in a broad range from 0 to 8−12 bonds, and the distribution becomes more symmetric with decreasing randomness of the sequence of reactive sites.
Because the effective potential V eff (r) is the free energy cost of changing the mutual distance from infinity to r, the difference between the effective potentials without and with intermolecular bonding is ΔβV eff (r) = βV eff,only intra (r) − βV eff,all bonds (r) = βΔU(r) − k B −1 ΔS(r), with ΔU(r) and ΔS(r) the corresponding energetic and entropic changes. For the same pair of polymers, switching intermolecular bonding on or off should not change excluded volume interactions significantly, and as mentioned before, it does not affect the total number of bonds. Therefore, ΔU(r) ≈ 0, and the difference between the effective potentials without and with intermolecular bonding is essentially of entropic origin, i.e., ΔβV eff (r) ≈ −k B −1 ΔS(r). Figure 6 shows (lines) the corresponding data for ΔβV eff (r). Because this quantity is positive for all distances, it is clear that forming intermolecular bonds involves an entropic gain with respect to the only intramolecularly bonded system. In principle, intermolecular bonds limit conformational and translational fluctuations, leading to an entropic loss. Therefore, there should be a source of entropic gain that exceeds the former loss, resulting in a net entropic gain when intermolecular bonds are formed. As can be seen in Figure 6, the net entropic gain is qualitatively given by k B times the number of intermolecular bonds; i.e., the number of additional states of the pair of polymers that are introduced by intermolecular bonding is essentially the exponential of the number of intermolecular bonds.
The mechanism leading to the observed entropic gain is not clear. The concept of combinatorial entropy, 44 accounting for the different connectivities of the bonding network, has been invoked to accurately describe a similar effect in the case of hard nanoparticles grafted by chains with sticky ends. An expression has been proposed for the number of bonding patterns that can be produced by the sticky ends that can, at each distance, potentially bind to the other nanoparticle. 45 Though it is plausible that the combinatorial entropy is a major contribution to the ΔβV eff shown in Figure 6, obtaining an Macromolecules pubs.acs.org/Macromolecules Article analytical accurate expression for our system is highly nontrivial 46 and is beyond the scope of this work.
On the other hand, as can be seen by comparing Figures 5  and 6, increasing the number of intermolecular bonds leads to lowering the effective potential. This is consistent with the fact that V eff (r)/k B T = −ln g(r), with g(r) the radial distribution of the centers-of-mass. 26 A higher number of intermolecular bonds leads to more tightly linked pairs, resulting in higher values of g(r) at short distances and, through the negative dependence, to lower values of V eff (r)/k B T. A similar trend should be found by increasing the total number of bonds (and concomitantly the intermolecular ones) through rising the ratio of the bond to the thermal energy.

CROWDED SOLUTIONS AND PHASE BEHAVIOR
The main motivation behind the coarse-graining approach is to reduce as much as possible the degrees of freedom that define the system. Deriving the expression of an effective potential V eff able to mimic the interactions between macromolecules enables the description of them only in terms of a few coordinates (usually the centers-of-mass). Thus, in a dense system as a crowded solution the degrees of freedom associated with the individual monomers are wiped out and the whole solution is effectively described as a fluid of particles interacting through the obtained V eff . This strategy largely reduces the computational cost of the all-monomer simulations, allows to investigate longer time and length scales, and facilitates the applications of methods from e.g. liquid state theory. However, it involves a strong assumption; namely, because the effective potential has been derived for two polymers in the absence of others, its use implicitly neglects the effective many-body interactions in the crowded solution. In general, this approximation is justified and works well for densities below the overlap concentration, but it fails, even severely, as one goes deep in the semidilute and concentrated regimes. 26 A well-known effect of the many-body interactions in dense solutions is the shrinkage found in simple linear chains, leading to the change from self-avoiding to Gaussian chain statistics. 21 Recent simulations of solutions of reversibly cross-linking linear chains similar to those investigated here have shown, interestingly, that the polymer size and shape are weakly affected by the concentration, essentially retaining the conformational properties of high dilution. 20 Instead of shrinking, the chains keep such mean conformations through forming a few intermolecular bonds with their neighbors. This weak effect of the concentration on the molecular conformations suggests that the many-body interactions experienced by a tagged couple of chains are in a first approximation given by a flat energy landscape. In such conditions the effective potential derived at high dilution may provide a good description of the structural properties of the solution even far above the overlap concentration. Figure 7 shows the radii of gyration, normalized by their values at ρ = 0, as a function of the normalized concentration ρ/ρ* for the linear chains and rings with reversible bonds, both in the pure systems and in the linear/ring mixture. The results for the pure linear case confirm those of the model of ref 20, with a shrinkage of just 4% at about 7 times the overlap concentration. A much steeper dependence on the concentration is found for the case of rings, with a shrinkage of 15% at the highest simulated concentration of about 4 times the overlap concentration (for comparison, at the same effective density the shrinkage of the linear chains is <2%). This very different response of the molecular size of linear chains and rings to crowding is also found in the mixture of both molecules, though differences are less pronounced than in the pure systems. In the mixture the size of the linear chains shows a steeper dependence on the concentration than in the pure system, whereas the rings show the opposite effect. Having said this, in all cases the shrinkage is much weaker that in the absence of bonding. Solid symbols in Figure 7 are the values for the unbonded precursors at the highest effective densities of the bonded counterparts. Such values have been estimated through the power laws R g /R g0 ∼ (ρ/ρ*) −1/8 (linear chains 21 ) and R g /R g0 ∼ (ρ/ρ*) −1/4 (unentangled rings 47 ). Shrinkage factors of 22% (linear) and 30% (ring) vs the respective aforementioned values of 4% and 15% are obtained, demonstrating the dramatic effect of intermolecular bonding on reducing the impact of crowding on the molecular conformations.
Beyond the effect of the concentration on the molecular size, the scattering form factor provides more detailed information about the molecular conformations. The form factor is calculated as where r jk = |r i − r j |; the sum is performed over all pairs of monomers j, k belonging to the same polymer and is averaged over all the polymers in the solution and different configurations. Figure 8 shows for both linear and ring architectures the form factor at high dilution and at the highest investigated concentration. As q grows, the form factor shows the crossover from the limit w(q = 0) = N m to the fractal regime, 21 w(q) ∼ q −1/ν , which originates from the scaling of the intramolecular distances with the contour length. Similar to the overall molecular size, we find a tiny effect of crowding on the effective exponent of the linear chains, which changes from ν = 0.58 to 0.54 from high dilution to ρ/ρ* ∼ 7, i.e., a narrow range between the Flory value (ν = 0.59) for self-avoiding chains and ν = 1/2 for Gaussian chains. The more pronounced effect of crowding on the molecular size of rings is also reflected in the scaling behavior, with a change from ν = 0.44 at high dilution to ν = 0.35 at ρ/ρ* ∼ 4, resembling crumpled globule behavior 22,23 (ν = 1/3). Similar trends are found in the 50/50 mixture of linear chains and rings. Consistently with the observations for the molecular size, the conformations of the linear chains in the mixture are slightly more affected by crowding than in the pure system, and the opposite effect is found for the rings.
In summary, Figures 7 and 8 show that the typical conformations of the linear chains are weakly distorted by crowding, and hence the two-body approximation under which the effective potential is derived might work reasonably even at unusually high densities, far above the overlap concentration. Comparatively, crowding has a stronger effect on the conformations of the rings, and their effective potential is expected to work in a narrower range of concentrations than in their linear counterparts. In what follows we test these expectations by comparing the results for the all-monomer solutions with those for the corresponding effective fluids. Moreover, we test the validity of mean-field DFT in our systems through calculations from test particle route (TPR). As mentioned in section 2, all the simulated solutions correspond to sequences of type "gap 1". The interactions in the all-monomer simulations are given by eqs 1−4. The data for the corresponding effective potentials of Figure 5 were fitted by the functions of Table 1, and these functions were used in the simulations and TPR calculations of the effective fluids. Namely, the "linear−linear, all bonds" and the "ring− ring, all bonds" potentials were used in the effective fluids of the pure (one-component) systems of linear and ring polymers with reversible bonds. They were also used for the linear− linear and ring−ring interactions in the effective linear−ring mixture, while the "linear−ring, all bonds" potential was used for the linear−ring interactions. In the mixture (A/B) of linear chains with orthogonal chemistry, the "linear−linear, all bonds" potential was used for the A−A and B−B interactions. Because by construction there were no intermolecular A−B bonds in the all-monomer simulations, the "linear−linear, only intra" potential was used for the A−B interactions in the effective fluid. Figure 9 shows the radial distribution function g(r) of the centers-of-mass in the pure solutions of linear chains with reversible bonds. Figure 9a compares the correlations for the centers-of-mass of the real all-monomer (AM) system with those for the particles of the effective fluid (EF). Figure 9b  compares the results for the effective fluid with the calculations from TPR. An excellent agreement between effective fluid and TPR is obtained, demonstrating the validity of the mean-field approximation for the effective fluid even at low densities. The comparison between the all-monomer and effective fluid reveals some interesting trends. Contrary to the usual observations in macromolecular systems, the effective potential provides a very good description of the real system at ρ/ρ* > 5, i.e., far above the overlap concentration, where many-body effects are usually expected. This finding confirms that the many-body effects are basically averaged out and lead to a flat energy landscape. Again contrary to the usual observations, there are systematic differences between the all-monomer and effective fluid at densities below the overlap concentration, even at values as low as ρ/ρ* ∼ 0.1, for which one might expect an excellent accuracy of the two-body approximation.  As can be seen in Figure 9a, the g(r) for the all-monomer system is shifted to longer distances, indicating less interpenetration than predicted by the effective fluid. The reason for this small but significant disagreement is likely the significant number of clusters of three polymers found at low concentrations in the real system. Figure S10 shows the cluster size distribution P(n) at the lowest investigated concentration, where n is the number of polymers in a cluster and two polymers belong to a same cluster if they are mutually linked by at least one intermolecular bond. As can be seen, the ratio of clusters of n = 3 vs those of n = 2 is non-negligible (about 0.1). In these clusters (which do not exist in simple systems with no bonds) the three-body interaction cannot be oversimplified by a flat landscape, and the two-body approximation just gives a semiquantitative description of the static correlations.
Results for the solutions of rings with reversible bonds are shown in Figure 10. In comparison with the linear case, the allmonomer rings show a larger correlation hole and therefore a weaker interpenetration. This is consistent with the observed stronger response of their conformations to crowding ( Figures  7 and 8), which leads to objects similar to crumpled globules (ν ∼ 1/3) and therefore less penetrable than their linear counterparts (ν ∼ 0.5). Although at low concentrations there is still a systematic small disagreement between the g(r) of the effective fluid and the all-monomer system, this effect is weaker than for the linear counterparts. This is consistent with the smaller number of three-body clusters found for the rings ( Figure S10). In this case the ratio of n = 3 vs n = 2 clusters is about 0.05. For concentrations higher than ρ/ρ* the effective fluid provides a much worse description than in the linear system, and indeed the all-monomer solution of rings does not show the peak at r = 0 found in the effective fluid. In a similar fashion to the simple case of rings without bonds, the peak formed at r = 0 and growing with the concentration is the signature of a fluid of clusters formed by strongly interpenetrated particles. The effective fluid will ultimately show a transition to a cluster crystal phase, where the clusters are arranged in the nodes of a regular lattice that is sustained through incessant hopping of the particles between the clusters. The existence of cluster crystal phases is predicted within mean-field DFT for potentials that are bounded and show negative values in their Fourier transform. 48 Both conditions are fulfilled by the effective potentials of the rings with reversible bonds. Indeed, they can be described by generalized exponential functions (Table 1), which for exponents higher than 2 have negative Fourier components. 33 Moreover, the mean-field approximation is justified, as can be seen in Figure 10b by the good agreement between the TPR and the simulations of the effective fluid. However, the cluster fluid is not found in the all-monomer system. As found for simple rings without bonding interactions, 35 the preferred crumpled globular conformations prevent the degree of nesting and threading needed to form the characteristic peak at r = 0. Figure 11 shows the partial correlations of the radial distribution function for the 50/50 ring−linear mixture. As in Figures 9 and 10 the left column compares AM and EF simulations, whereas the right column compares the EF simulations with the predictions of TPR. The top panels (a, b) display the partial correlations between the linear chains, g ll (r). The middle panels (c, d) show the cross-correlations (linear− ring, g lr (r)), and the correlations between the rings (g rr (r)) are displayed in the bottom panels (e, f). At low and moderate concentrations, the small but systematic deviations between the AM and EF for the linear−linear correlations in the mixture are similar to those found in the pure system. However, whereas at large concentrations (ρ/ρ* > 4) there is a very good agreement between the AM simulations and EF in the pure linear system, significant differences are observed in the mixture. This suggests that the picture of an effective flat energy landscape describing the many-body interactions in the pure linear system is an oversimplification when the linear neighbors are partially substituted by rings adopting less penetrable crumpled globular conformations and hence leading to an heterogeneous landscape. This is consistent with the found deviations between the EF and TPR (see panel (b)), in contrast to the excellent agreement observed in the pure system. The description of the ring−ring AM correlations by the effective potentials is improved in the mixture with respect to the pure solutions. Indeed, the presence of a 50% of particles (representing the linear chains) in the EF interacting through Gaussian potentials (which do not lead to cluster phases) reduces the tendency to the cluster phase of the particles representing the rings, and the EF becomes closer to the real system where no peak at r = 0 is found. On the other hand, the TPR provides a worse description of the ring−ring correlations in the EF of the mixture than in the EF of the pure system of rings. Again, this might be related to the structural heterogeneity of the EF of the mixture that worsens the meanfield approach of TPR, though surprisingly, TPR does provide a very good description of the cross-correlations (linear−ring) in the EF. A reasonably good agreement is also found between the cross-correlations in the AM and EF systems. The results for the cross-correlations in panel (c) show a good mixing of the linear chains and rings with reversible bonds, with no signatures of segregation or incoming phase separation. Indeed, the correlation holes are just intermediate between those for the self-correlations.
Finally, we push further our investigation on the validity of effective potentials to describe correlations in crowded solutions of polymers with reversible bonds. Our last question is whether it is possible in our model to form interpenetrated networks (IPNs) from two polymers with reversible bonds and orthogonal chemistry. For this purpose, we consider a linear binary mixture with the same fraction x = 50% for both components. As mentioned before, in the all-monomer simulations the WCA, FENE, and reversible bonding interactions are identical for both components, with the only difference that intermolecular bonding is switched off between chains of different components. In the effective fluid, the interactions between particles of the same component are the same as in the EF of the pure linear case ("linear−linear, all bonds" in Table 1), whereas for the cross-interactions we use the effective potentials derived in the absence of intermolecular bonding ("linear−linear, only intra"). To have a first idea of the emerging scenario for this system, we obtain the theoretical phase diagram for the effective fluid in the plane of reduced concentration (ρ/ρ*) vs composition (0 ≤ x ≤ 1) of the mixture using the random phase approximation for the partial correlations as a closure to the Ornstein−Zernike relation. 29,42,49,50 This should be a reasonable approximation on the basis of the observed quality of the mean-field TPR. We find a spinodal line (dashed line in Figure 12) attesting to the existence of a region with macrophase separation (demixing) which, because self-interactions for both components are identical, becomes symmetric with respect to the composition. On the other hand, we find that, except for very asymmetric compositions, the system demixes when the density is increased slightly above the overlap concentration. Because the onset of network percolation occurs at such concentrations or above them, 20 the theoretical phase diagram of Figure 12 suggests that it is not possible to form an IPN in the mixture of chains with reversible bonds, this being frustrated by the demixing of both components. Figure 13 shows AM and EF simulation snapshots at different concentrations for the 50/50 mixture of linear chains with reversible bonds and orthogonal chemistry. The beads represent the monomers and the effective particles in the AM and EF systems, respectively. By depiction of the two components with different colors, demixing (which as anticipated in section 2 occurs spontaneously by evolution from an initially mixed state) is evident and confirms the expectation from the theoretical phase diagram. This is quantitatively reflected in the partial components of the total static structure factor of the molecular centers-of-mass, S αβ (q), where α, β refer to the components (1,2) of the mixture, so that S 11 (q) and S 22 (q) represent correlations within a same component and S 12 (q) represents cross-correlations between chains of different components. These quantities are calculated as In this equation N α is the number of relevant coordinates of the α-component in the simulation box (the molecular centersof-mass in the AM and all the effective particles in the EF), and r j α denotes the coordinate of the jth molecule of the αcomponent. The average is performed over several realizations of the box and different runs at the same concentration. The total structure factor, S(q), accounting for all the correlations without distinguishing components of the mixture, is just obtained by running the sum over all pairs of coordinates in the box (irrespective of their respective components) and normalizing the sum by the inverse of the total of number molecules N A + N B . Figure 14 shows the total S(q) (panel (a)) and the partial components S αβ (q) (panels (b−d)) of the molecular centers-of-mass in the AM system. It should be noted that because the composition is equimolar and the selfinteractions of the two components are identical, S 11 (q) = S 22 (q). No signatures of growing length scales are observed in the total S(q), which shows the typical behavior of a homogeneous fluid with increasing the concentration. The growing length scales of the two separating phases are evidenced by the growing peaks of the partial S 11 (q) =   Phase separation has been found in a simplified mixture where all nonbonded and bonded interactions are identical, with the only constraint that intermolecular bonds between different components of the mixture are not allowed. A more realistic model should at least introduce different activation energies for the two kinds of orthogonal reactive sites. This would likely break the symmetry of the phase diagram with respect to the mixture composition, but the qualitative emerging scenario (demixing and impossibility of forming the IPN in equilibrium) is robust. Indeed, demixing is inherently connected to the more repulsive character of the cross-interactions than of the self-interactions. Introducing different bond activation energies will lead to different bonding probabilities and hence different self-interactions of the two components, but the cross-interaction should still be much more repulsive than the self-interactions because, as discussed in section 3, the latter contain the combinatorial entropic gain associated with intermolecular bonding, this being absent between polymers of different components of the mixture. Having said this, there is plenty of evidence in the literature on formation of IPNs with purely reversible cross-links. 51−55 Our results suggests that such IPNs are kinetically trapped states. The typical activation energies of the dynamic bonds in these IPNs are of several hundreds of kJ/mol, 55 i.e., of the order of 100k B T, whereas in our simulations they are about 10k B T. Moreover the experimental chains are much longer than the unentangled chains used here (the entanglement monomer density for our linear precursors is 20 ρ e ≳ 0.42σ −3 ). Thus, our results suggest that in real systems the combination of both high molecular weights and long lifetimes of the bonds creates large barriers impeding relaxation to the equilibrium demixed state, and the IPNs (created in out-of-equilibrium conditions) remain stable.

CONCLUSIONS
We have systematically investigated effective potentials between polymeric molecules functionalized with groups that can form intra-and intermolecular reversible bonds. A rich scenario emerges for the dependence of the effective potential on the relevant control parameters. In spite of the additional complexity introduced by the high number of instantaneous intramolecular loops originated by the reversible cross-links, the topological interaction of the unbonded precursor (linear or ring) still has a dominant contribution in the bonded state, leading to very different strengths of the effective interaction (being more repulsive for the ring-based system). Even if the molecular weight and the fraction of reactive sites are fixed, the effective potentials exhibit a significant dependence on the degree of randomness of the sequence of reactive sites (from fully random to periodic). If the reactive sites of the two polymers are orthogonal, so that only intramolecular bonds are formed, decreasing randomness leads to longer intramolecular loops, which hinders interpenetrability and leads to a stronger effective intermolecular repulsion. The opposite effect is found if reactive sites of both polymers are identical and both intraand intermolecular bonding occur. We suggest that the free energy loss caused by the intermolecular bonds is mainly given by combinatorial entropy arising from the exponential number of bonding patters that the two intermolecularly bonded polymers can adopt.
We have explored the accuracy of the effective potentials to describe the equilibrium correlations between centers-of-mass in the crowded solutions. In the case of the linear chains a very good agreement between the effective fluid and the allmonomer simulations is found ever far above the overlap concentration. This is consistent with the fact that shrinking is highly prevented by forming intermolecular bonds with neighboring chains, which makes the conformations at high dilution weakly sensitive to crowding, and many-body effects basically contribute as a flat energy landscape. In a similar fashion to the case of rings with no bonds, the comparison with the effective fluid is less satisfactory in the system of rings with reactive sites, which does not show the cluster phase predicted by the effective fluid. This is consistent with the crowdingdriven collapse to crumpled globule-like conformations, reflecting the relevance of the many-body interactions. We have further extended our investigation to a 50/50 mixture of the former types of polymers. The results for the partial correlations are qualitatively similar to those of the pure polymers and the system is fully miscible.
Finally, we have explored the possibility of forming two interpenetrated networks in a linear−linear mixture where the reactive sites of the two components are orthogonal; i.e., intermolecular bonds only occur between chains of the same component. In agreement with the energetic penalty found for the effective cross-interaction potential, the simulations of the effective fluid, and the phase diagram obtained by the test particle route, no interpenetrated networks are found, and the two components demix. This result suggests that real interpenetrated networks, where the lifetimes of the reversible bonds are much longer than in our simulations, are kinetically trapped states with large entropic barriers impeding the relaxation to the equilibrium demixed state (arrested demixing). Our results may motivate future experimental tests in mixtures of oligomers with low bond energies. On the other hand, an interesting problem to address in the future is the accuracy of the effective fluid approach in dual networks, where both types of orthogonal reactive sites are present in all the chains, including the determination of the phase behavior in mixtures with different fractions of both sites in each component. Another future line of research would be to improve the description of the real system through the incorporation of additional degrees of freedom. A promising approach 56 is to introduce an intermediate pair potential that depends on the instantaneous values of the intermolecular distance and of the two molecular sizes. The latter are coupled to the density of the solution through the equations of motion, leading to an effective potential (averaged over the size distribution) that becomes density dependent. Work in these directions is in progress.
Derivation of TPR equations; Tables S1 and S2: details of the simulated boxes for the crowded solutions; Figures S1 and S2: sequence of reactive sites for the systems used in the derivation of the effective potentials; Figures S3−S10: characterization of several conformational and structural properties (PDF)